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Abstract This paper is a demonstration of liow tiie WKB approximation can be 
used to help solve the linearised 3D MHD equations. Using Charpit's Method and a 
Runge-Kutta numerical scheme, we have demonstrated this technique for a potential 
3D magnetic null point, B = {x, ey — {e + 1) z). Under our cold plasma assumption, 
we have considered two types of wave propagation: fast magnetoacoustic and Alfven 
waves. We find that the fast magnetoacoustic wave experiences refraction towards the 
magnetic null point, and that the effect of this refraction depends upon the Alfven 
speed profile. The wave, and thus the wave energy, accumulates at the null point. 
We have found that current build up is exponential and the exponent is dependent 
upon e. Thus, for the fast wave there is preferential heating at the null point. For 
the Alfven wave, we find that the wave propagates along the fieldlines. For an Alfven 
wave generated along the fan-plane, the wave accumulates along the spine. For an 
Alfven wave generated across the spine, the value of e determines where the wave 
accumulation will occur: fan-plane (e = 1), along the x— axis (0 < e < 1) or along 
the y— axis (e > 1). We have shown analytically that currents build up exponentially, 
leading to preferential heating in these areas. The work described here highlights the 
importance of understanding the magnetic topology of the coronal magnetic field for 
the location of wave heating. 

Keywords: Magnetohydrodynamics; Waves, Propagation; Magnetic fields, Models; 
Heating, Coronal 

1. Introduction 

The WKB approximation is an asymptotic approximation technique which can be 
used when a system contains a large parameter (see e.g. Bender and Orszag, 1978). 
Hence, the WKB method can be used in a system where a wave propagates through a 
background medium which varies on some spatial scale which is much longer than the 
wavelength of the wave. The SOHO and TRACE satellites have recently observed 
MHD wave motions in the corona, i.e. fast and slow magnetoacoustic waves and 
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Alfven waves (see reviews by Nakariakov and Verwichte, 2005; De Moortel, 2005; 
2006) . The coronal magnetic field plays a fundamental role in their propagation and 
to begin to understand this inhomogeneous magnetised environment, it is useful to 
look at the structure (topology) of the magnetic field itself. Potential-field extrapo- 
lations of the coronal magnetic field can be made from photospheric magnetograms. 
Such extrapolations show the existence of an important feature of the topology: null 
points. Null points are points in the field where the magnetic field, and hence the 
Alfven speed, is zero. Detailed investigations of the coronal magnetic field, using such 
potential field calculations, can be found in Beveridge, Priest, and Brown (2002) and 
Brown and Priest (2001). 

McLaughlin and Hood (2004) found that for a single 2D null point, the fast 
magnetoacoustic wave was attracted to the null and the wave energy accumulated 
there. In addition, they found that the Alfven wave energy accumulated along the 
separatrices of the topology. They solved the 2D linearised MHD equations numer- 
ically and compared the results with a WKB approximation: the agreement was 
excellent. From their work and other examples {e.g. Galsgaard, Priest, and Titov, 
2003; McLaughhn and Hood, 2005; 2006a; Khomenko and CoUados, 2006) it has 
been clearly demonstrated that the WKB approximation can provide a vital link 
between analytical and numerical work, and often provides the critical insight to 
understanding the physical results. This paper demonstrates the methodology of 
how to apply the WKB approximation in linear 3D MHD. We believe that with 
the vast amount of 3D modelling currently being undertaken, applying this WKB 
technique to 3D will be very useful and beneficial to modellers in the near future. 

The work undertaken by Galsgaard, Priest, and Titov (2003) deserves special 
mention here. They performed numerical experiments on the efi^ect of twisting the 
spine of a 3D null point, and described the resultant wave propagation towards 
the null. They found that when the fieldlines around the spine are perturbed in a 
rotationally symmetric manner, a twist wave (essentially an Alfven wave) propagates 
towards the null along the fieldlines. Whilst this Alfven wave spreads out as the null 
is approached, a fast-mode wave focuses on the null and wraps around it. They 
concluded that the driving of the fast wave was likely to come from a non-linear 
coupling to the Alfven wave (Nakariakov, Roberts, and Murawski, 1997). They also 
compare their results with a WKB approximation and find that, for the (3 = fast 
wave, the wavefront wraps around the null point as it contracts towards it. They 
perform their WKB approximation in cylindrical polar coordinates and thus their 
resultant equations are two-dimensional (since a simple 3D null point is essentially 
2D in cylindrical coordinates). In contrast, we solve the WKB equations for three 
Cartesian components, and thus we can solve for more general disturbances and more 
general boundary conditions. This also allows us to concentrate on the transient 
features that are not always apparent when only cylindrically symmetric solutions 
are permitted. 

More recently, Pontin and Galsgaard (2007) and Pontin, Bhattacharjee, and Gals- 
gaard (2007) have performed numerical simulations in which the spine and fan of 
a 3D null point are subject to rotational and shear perturbations. They found that 
rotations of the fan plane lead to current sheets in the location of the spine and 
rotations about the spine lead to current sheets in the fan. In addition, shearing 
perturbations lead to 3D localised current sheets focused at the null point itself. 
This general behaviour is in good agreement with the work presented in this paper. 
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i.e. current accumulation at certain parts of the topology. However, the primary mo- 
tivation in Pontin and Galsgaard (2007) and Pontin, Bhattacharjee, and Galsgaard 
(2007) was to investigate current-sheet formation and reconnection rates, whereas 
the techniques described in this paper focus on MHD wave-mode propagation and 
interpretation. 

The propagation of fast magnetoacoustic waves in an inhomogeneous coronal 
plasma has been investigated by Nakariakov and Roberts (1995), who showed how 
the waves are refracted into regions of low Alfven speed. In the case of null points, 
the Alfven speed actually drops to zero. 

The paper has the following outline: In Section 2, the basic equations are de- 
scribed. Section 3 details the 3D WKB approximation utilised in this paper. The 
results for the fast wave and Alfven waves are shown in Section 4 and 5. The 
conclusions and discussion are presented in Section 7. There are four appendices 
which complement the results in the main text. 



2. Basic Equations 

The usual resistive, adiabatic MHD equations for a plasma in the solar corona are 
used: 



dB 

'dt 



V X (v X B) + r7V^B , (2) 



If +V-(pv) = , (3) 

^ + V • Vp = -7pV.v , (4) 
= V X B , (5) 

where v is the plasma velocity, p is the mass density, p is the gas pressure, B is 
the magnetic induction (usually called the magnetic field), j is the electric current, 
g is gravitational acceleration, 7 is the ratio of specified heats, 77 is the magnetic 
diffusivity and fi is the magnetic permeability. 



2.1. Basic Equilibrium 



We choose a 3D magnetic null point for our equilibrium field, of the form: 

Bo = j{x,ey,-{e + l)z) , (6) 

where B is a characteristic field strength, L is the length scale for magnetic field 
variations and the parameter e is related to the predominate direction of alignment 
of the fieldlines in the fan plane. Parnell et al. (1996) investigated and classified the 
different types of linear magnetic null points that can exist (our e parameter is called 
p in their work). Topologically, this 3D null consists of two key parts: the 2;— axis 
represents a special, isolated fieldline called the spine which approaches the null from 
above and below (Priest and Titov, 1996) and the xy— plane through z = is known 
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Figure 1. Left: Proper radial null point, described by B = {x,y, —2z), i.e. e = 1. Right: 
Improper radial null point, described by B = (x, ^y,—^z), i.e. e = ^. Note for e = ^, the 
field lines rapidly curve such that they run parallel to the x— axis along y = 0. In both figures, 
the z— axis indicates the spine and the xy— plane at 2: = denotes the fan. The red fieldlines 
have been tracked from the 2 = 1 plane, the blue from z = — 1. 



as the fan and consists of a surface of fieldlines spreading out radially from the null. 
Figure 1 shows two examples of 3D null points: e = 1 (left) and e = 1/2 (right). 
Titov and Hornig (2000) have investigated the steady state structures of magnetic 
null points. 

Equation (6) is the general expression for the linear field about a potential mag- 
netic null point (Parnell et al, 1996: Section IV). In this paper, we only consider 
e > and so all nulls we describe are positive nulls, i.e. the spine points into the null 
and the field lines in the fan are directed away. In addition, all potential nulls are 
designated radial, i.e. there is no spiral motions in the fan-plane. In general, there 
are three cases to consider: 

• e = 1: describes a proper null (Figure 1: Left). This magnetic null has cylindrical 
symmetry about the spine axis. 

• e > 0, e ^ 1: describes an improper null (Figure 1: Right). Field lines rapidly 
curve such that they run parallel to the a;— axis if < e < 1 and parallel to the 
y— axis if e > 1. 

• e = 0: equation (6) reduces to the X-point potential field in the xz— plane and 
forms a null line along the y— axis through x = 2 = 0. MHD wave propagation 
in this 2D configuration has been studied extensively by McLaughlin and Hood 
(2004; 2005; 2006a). 



2.2. Assumptions and Simplifications 

In this paper, the linearised MHD equations are used to study the nature of wave 
propagation near the null point. Using subscripts of for equilibrium quantities and 
1 for perturbed cjuantities. Equations (l)-(5) become: 

= -Vpi +jo X Bi -Fji X Bo + pig , (7) 
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(10) 
(11) 



(8) 



(9) 



We now consider several simplifications to our system. We will only be considering 
a potential equilibrium magnetic field (V x Bo = 0) in an ideal system (rj = 0). We 
will also assume the equilibrium gas density (po) is uniform. A spatial variation in po 
can cause phase mixing (Heyvaerts and Priest, 1983; De Moortel et al, 1999; Hood, 
Brooks, and Wright, 2002). In addition, we ignore the effect of gravity on the system 
{i.e. we set g = 0). Finally, in this paper we assume a cold plasma, i.e. Cs = 

\/lPo/po = 0. 

We will not discuss Equation (9) further as it can be solved once we know vi. In 
fact, under tlie assumptions of linearisation and no gravity, it has no influence on the 
momentum equation and so in effect the plasma is arbitrarily compressible (Craig 
and Watson, 1992). 

We now non-dimensionalise the above equations as follows: let vi = vv^j^, Bq = 
-BBq, Bi = BBl, X = Lx* , z = Lz* , V = -^V* and t = tt* , where we let * denote 
a dimensionless quantity and v, B, L, and t are constants with the dimensions of 
the variable that they are scaling. In addition, po and po are constants as these 
equilibrimn quantities are uniform (i.e. Po = Po = !)• We then set B/y/fipo = v and 
V = L/t (setting v as a constant background Alfven speed). Under these scalings, 
t* = I (for example) refers to t = t = L / v; i.e. the (background) Alfven time taken 
to travel a distance L. For the rest of tliis paper, we drop the star indices; tlie fact 
that they are now non-dimensionalised is understood. 

These non-dimensionalised equations can be combined to form one single equa- 
tion: 



3. WKB Approximation 

In this paper, we will be looking for WKB solutions of the form: 



where a is a constant. In addition, we define a; = 0t as the frequency and k = V^A = 
{<px, <py, 0z) = (p, q, r) as the wavevector. (/), and its derivatives, are considered to be 
the large parameters in our system. 

One of the difficulties associated with 3D MHD wave propagation is distinguishing 
between the three different wave types, i.e. between the fast and slow magnetoacous- 
tic waves and the Alfven wave. To aid us in our interpretation, we now define a new 
coordinate system: (Bq, k, Bq x k), where k is our wavevector as defined above. This 
coordinate system fully describes all three directions in space when Bq and k are not 
parallel to each other, i.e. k ^ ABq, where A is some constant of proportionality. In 



—VI = {V X [V X (vi X Bo)]} X Bq 



(12) 
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the work below, we will proceed assuming k ^ ABq. The scenario where k = ABq is 
looked at in Appendix A. In fact, the work described below is also valid for k = ABo 
with the consequence that the solution is degenerate, i.e. the waves recovered are 
identical and cannot be distinguished. 

We now substitute v = ae"^^"^'^'^'*^ into Equation (12) and make the WKB 
approximation such that cf> ^ 1. Taking the dot product with Bq, k and Bo x k 
gives three velocity components which, in matrix form, are: 

uj"^ 1/v-Bo\/o\ 

(Bo-k)|k[2 L^2_|B^|2|j^|2 Q ^.^ = 

- (Bo • k)^ J V V ■ Bo X k / \0 J 

The matrix of these three coupled Equations must have zero determinant so as not 
to have a trivial solution. Thus, taking the determinant gives: 

T{(l),uj,t,Bo,k) = (c^'-o) - |Bo|' Ikj^) (^cj^ - {Bo ■ kf) =0 , (14) 

where .F is a first-order, non-linear PDE. Equation (14) has two solutions, corre- 
sponding to two different MHD wave types (in general three, but the slow wave has 
vanished under the cold plasma approximation). The two solutions correspond to 
the fast magnetoacoustic wave and to the Alfven wave. 

In Sections 4 and 5, we will examine each of these wave solutions in detail for the 
3D magnetic null point configuration described by Equation (6) for both e = 1/2 
and e = 1. However, it should be noted that the technique described above is valid 
for any 3D magnetic configuration. The case where the two roots of Equation (14) 
are the same is examined in Appendix A. 



4. Fast Wave 

Let us first consider the fast wave solution, and hence we assume 7^ (Bo • k)^. 
Thus, Ecjuation (14) simplifies to: 

= - (x" + + (e + if z") (j? + g2 + = (15) 

We can now use Charpit's Method (see e.g. Evans, Blackledge, and Yardley, 1999) to 
solve this first-order PDE, where we assume our variables depend upon some inde- 
pendent parameter s in characteristic space. Charpit's Method replaces a first-order 
PDE with a set of characteristics that are a system of ODEs. Charpit's Equations 
take the form: 

^ = LA+k-— V — = —T ^ = A^ 

ds Y duj 9k J ' ds duj ' ds 9k ' 

duj _ _( d 9\ dk__/9 .d\ 

ds " ydt^'^dct))^ ' ds " \d^^ d4>) ' 

where, as previously defined, k = {p,q,r) and x = {x,y,z). In general, the coupled 
Equations have to be solved numerically, but analytical solutions have been found in 
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2D (McLaughlin and Hood, 2004). These ODEs are subject to the initial conditions 
(j) = (j)Q[s = 0), X = xo{s = 0), y = yo(s = 0), 2 = zq{s = 0), t = to{s = 0), 
P = Po{s = 0), q = qo{s = 0), r = ro{s = 0), and u) = ujo{s = 0) and, in the 
following work, are solved numerically using a fourth-order Runge-Kutta method. 

In addition, note that there are no boundary conditions in the usual sense: the 
variables are solved using Charpit's Method (essentially tire method of character- 
istics) and the resulting characteristics are only dependent upon initial position 
{xo,yo, zo,to) and distance travelled along the characteristic; s. Thus, there are 
no computational boundaries and no boundary conditions (only initial conditions). 
In tliis paper, we have chosen to illustrate our results in tire domain — 1 < a: < 1, 
— I<j/<1,— 1<2:<1, and this choice is arbitrary. That the WKB solutions 
are independent of boundary conditions is actually an advantage over traditional 
numerical simulations; where the choice of boundary conditions can play a significant 
role. 

For the fast wave (Equation 15) Charpit's Equations are: 

dd> dt dx . dy . dz 

-f = , ^ = uj , — = -pA, -f = -qA, :r = 

ds ds ds ds ds 

— = , — = xB , — =e yB , — = [e + 1) zB (16 
ds ds ds ds 

where A = x^ + e^y^ + (e + 1)^ and B = + q^ + r^. 

From these Equations, we note tliat (j) = constant = <j)o and lu = constant = uq, 
i.e. constant frequency. In addition, t = los + to, where we arbitrarily set = 0, 
which correponds to the leading edge of the wave pulse starting at t = when s = 0. 
We can also construct the integral: 

— (xp + yq + zr) =0 xp -\- yq -\- zr = constant = xopo + 2/0?0 + ■ZO^'O (17) 

However, we are unable to find a second conserved quantity. 



4.1. Planar fast wave starting at zq = 1 



We now solve Equation (16) subject to the initial conditions: 

(/)o = 0, LUQ = 2-1: , -l<a;o<l, -l<yo<l, 20 = 1, 

po = , go = , ro = uJol^xl + e-^yl + {e + lf zl , (18) 

where we have arbitrarily chosen u>o = 27r and (^0 = 0. These initial conditions 
correspond to a planar fast wave being sent towards the null point from our upper 
boundary (along 2 = 20). 

Let us initially consider e = 1 (corresponding to Figure 1: Left). In Figure 2, we 
have plotted surfaces of constant which can be thought of as defining tire position 
of the wavefront, at various times. Since t = ujs, these correspond to different values 
of the parameter s, which quantifies distance travelled along the characteristic curve. 
We can clearly see that the fast wave experiences a refraction effect towards the null 
point, i.e. propagation towards regions of lower Alfven speed. A similar refraction 
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Figure 2. (e = 1) Surfaces of constant at four values of t, showing the behaviour of the 
(initially planar) wavefront that starts at —1 < x < 1, —1 < y < and z = 1 (blue) and 
— 1 <x<l, 0<y<l and z = I (red). The (arbitrary) colouring has been added to aid the 
reader in tracking the wave behaviour. This figure is also available as an mpg animation. 



effect was also seen in the 2D case (McLaughlin and Hood, 2004; 2006a). Thus, the 
fast wave is deformed from its initial planar profile. The wave, and hence all the 
wave energy, eventually accumulate at the null point. 



We can also use our WKB solution to plot the ray paths of individual fiuid 
elements from the initial wave. In Figure 3, we can see the ray paths for fluid elements 
that begin at four different starting points in the zq = 1 plane. These four ray paths 
are typical of the behaviour of the fast wave fiuid elements (more examples can be 
seen in the associated mpg movie). We can clearly see the refraction effect wrapping 
the fast wave fluid elements around the null point. 

Note that the magnitude of the refraction effect is the same for each fluid element 
that starts at the same radius from z = zq. Thus, the sy— plane projections are 
always straight lines. This is because the Alfven speed is the same for elements 
starting at the same radius from z = zo, i.e. VA{xo,yo, zo) = Xq + + Az'^, and 
the behaviour of the fast wave is entirely dominated by the Alfven speed profile. For 
e = 1, isosurfaces of Alfven speed form prolate spheroids (parallel to the spine). 

Let us now extend our study to improper null points. In Figure 4, we can see 
the ray paths for fluid elements that begin at the same four starting points as in 
Figure 3, but now for the magnetic field configuration seen in Figure 1: Right, i.e 
e = 1/2. The first thing to note is that the refraction effect still occurs in this 
configuration, as expected, and that the fiuid elements still eventually accumulate 
at the null point. However, the individual ray paths are different to those for e = 1. 
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Figure 3. (e = 1) Ray paths for fluid elements that begin at (xq, j/o,^o) = {^) (0, 1, 1), {B) 
(1, 0, 1), (C) (—1, 1, 1) and (D) (—0.5, 0.5, 1). This figure is also available as an mpg animation 
showing all —1 < xq < 1, zq = 1 along y = —1 and y = I, and all —1 < J/o !S li ^0 = 1 along 
X = —1 and X = I. Here, red indicates the 3D ray path and black indicates the xy, yz and xz 
projections of this ray path onto the respective planes. The blue dot indicates the position of 
the magnetic null point. 

Comparing panels 4A and 4S with 3 A and 3B, we see that along x = or y = 0, 
the ray paths are very similar. However, the distance travelled by the fluid element 
is actually different in the e = 1/2 simulation since the Alfven speeds, and hence 
the magnitude of the refraction, has changed. In fact, the refraction is weaker along 
y = (since for \B\ = + e^y^ + (e + 1)^ , \B\^^i > \B\^^i/2 ^long y = 0) and 
so the fluid element travels a longer distance than the equivalent e = 1 fluid element. 
Along a; = 0, the effect is more complicated, with l-Bl^^^ ^ l-^lc=i/2 only true for 

l^l > \y\- 

The differences in the ray paths are much more obvious when comparing panels 
AC and AD with 3C and 3Z). For e = 1/2, we see that the ray paths are now 
"corkscrew" spirals. This is because the Alfven speed profile is now varying in three 
directions, whereas for e = 1 the Alfven speed essentially varies in two directions: 
r = \J ^ and z. Thus, the xy— plane projections are no longer straight lines. 
For e = 1/2, isosurfaces of Alfven speed form scalene ellipsoids. 

4.2. Planar fast wave starting at yo = 1 
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Figure 4. (e = 1/2) Ray paths for fluid elements that begin at (xq, j/q, zq) = (A) (0, 1, 1), {B) 
(1, 0, 1), (C) (—1, 1, 1) and (D) (—0.5, 0.5, 1). This figure is also available as an mpg animation 
showing many more starting points. Red indicates the 3D ray path and black indicates the xy, 
yz and xz projections of this ray path onto the respective planes. The blue dot indicates the 
position of the magnetic null point. 

We again solve Equation (16) but now subject to tlie initial conditions: 
00 = , ti;o = 27r , -1 < a;o < 1 , yo = 1 , -1 < ^0 < 1 , 
po = , go = ujo/\lxl + e'^yl + {e+lf zl , tq = , (19) 

These initial conditions correspond to a fast wave being sent in from the side bound- 
ary (along y = yo). This choice of planar fast wave is incident perpendicular to the 
spine (z— axis). 

Let us first consider e = 1. Figure 5 shows surfaces of constant (f> at four values of t, 
showing the behaviour of the (initially planar) wavefront that starts at — 1 < a; < 1, 
y = I and —l<z<l. Again, we see the deformation of the (initially planar) 
wave due to the refraction effect and, again, the wave accumulates at the null point. 
However, the nature of this refraction is different to that seen in Figure 2, since the 
refraction varies in magnitude in different planes. We see the wavefront is initially 
"pinched" preferentially in the yz— plane (since Vj^{y,z) > Vj^{x,y) for \2z\ > |a::|). 

In Figure 6, we can see the ray paths for fluid elements that begins at four different 
starting points in the yo = 1 plane. Again, we can clearly see the refraction wrapping 
the fast wave elements around the null point, and the ray paths accumulate at the 
null point. 

We do not show the ray paths corresponding to a planar fast wave starting at 
a;o = 1 and approaching the null point, since these the ray paths behave identically to 
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Figure 5. (e = 1) Surfaces of constant at four values of t, showing the behaviour of the 
(initially planar) wavefront that starts at —1 < x < 1, y = 1 and —1 < z < (blue) and 
— 1 < X < I, y = 1 and < z < 1 (red). The (arbitrary) colouring has been added to aid the 
reader in tracking the wave behaviour. The yellow dot indicates the position of the magnetic 
null point. 

those in Figure 6 under the transformation (x, y) — > {—y, x) (for e = 1 configuration). 
The ray paths corresponding to a planar fast wave starting at yo = 1 and starting 
at xo = 1 in tire e = 1/2 magnetic configuration can be found in Figures 11 and 12 
in Appendix B. 

Thus, Sections 4.1 and 4.2 have shown that the fast wave experiences a refraction 
effect in the neighbourhood of a 3D magnetic null point and that in all of these 
cases, the main result is the same: the ray paths accumulate at the null point. Of 
course, the actual paths taken vary depending upon initial conditions and choice of 
e. Hence, we conclude that the fast wave, and thus the fast-wave energy, eventually 
accumulate at the 3D null point for all e and all initial conditions that generate a 
wave approaching the null. 

Finally, it should be noted that the behaviour of the fast wave is entirely domi- 
nated by the Alfven-speed profile, and since the magnetic field drops to zero at the 
null point, the wave will never actually reach there. However, there is still current 
accumulation and hence non-ideal effects may be able to extract the wave energy in 
a finite time. This is investigated in the next section. 

4.3. Current build up 

From Section 3, we know that for the fast wave v ■ Bo = and v • (Bq x k) = 0. 
This gives us two Equations for the three velocity variables: 

xvx + eyvy - (e + 1) zvz = , 
[eyr + (e + 1) zq] Vx - [xr + (e + 1) zp] Vy + [xq - eyp] Vz = ■ 

Thus, we can express two of the velocity components in terms of the third. 
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Figure 6. (e = 1) Ray paths for fluid elements that begin at {xq, yo, zo) = {^) (0, 1, 1), [B) 
(— 1, 1, 0), (C) (—1, 1, 1) and (D) (—0.5, 1, 0.5). This figure is also available as an mpg animation 
showing many more starting points. Red indicates the 3D ray path and black indicates the xy, 
yz and xz projections of this ray path onto the respective planes. The blue dot indicates the 
position of the magnetic null point. 



Recall from Section 2 that the perturbed electric current is given byji = VxBi. 
Thus, 

V X (v X Bo) —uj'Bi = k X (v X Bo) , 
ik X Bi = — ik X [k X (v X Bo)] /uj = i |k|^ (v x Bo) /lo , 




where we have made use of v • (Bq x k) = 0. From Equation (15), we can substitute 
for |k[ to obtain: 



Jl 



iJi 



iLu (v X Bo)/|Bor 

. [- (e + 1) zvy - eyvz, (e + 1) zvx + xvz,£yvx - xVy] 



rj.2 _|_ ^2yl + (g + \y ^2 



u: |v| / |Bo| 



(20) 



where we have used v-Bq = to simplify |v x Bo|. Thus, since |v| is bounded (from 
our assumed form of v seen in Equation 13) we can see that the current associated 
with the fast wave will grow ~ 1/ |Bo|. Equivalent behaviour was found for the fast 
wave in the 2D case (McLaughlin and Hood, 2004), i.e. e = 0. 
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Moreover, we can place limits on the magnitude of the current build up. From 
Equation (37) in Appendix C, we can place limits on |ji| such that: 

(e + l)i?0 -l-''l-ei?o ' ^ ^ 

where i?Q = Xg + j/q + 2:0 s-nd a = xopo + yoqo + zoro (see Appendix C) and we have 
assumed < e < 1. Thus, we can see that the current build up is bounded by two 
exponentially growing functions. 

We now demonstrate this current build up for two particular cases. Firstly, con- 
sider a planar fast wave starting at = 1 (Section 4.1). Here, we can solve Charpit's 
Equations for the fast wave (Equation 16) analytically for the initial conditions 
Pq = Qo = xq = yo = 0, i.e. along x = y = which is the path along which we 
expect the maximum current build up to occur. Under these conditions. Equation 
(16) reduces to: 

x = , y = , ^ = -{e + lfrz\ p = , q = , ^ = {e + if zr^ , 
as as 

where we have used initial conditions (18). We also note that our conserved quantity 
(Equation 17) states zr = zqvo = ujo / {e + 1). Thus: 

z = zoe-(^+i)"°^ = zoe-(^+^)* , r = roe^'+i)""^ = , ^\ e(^+^)* . (22) 

(e+ 1)^0 

As mentioned previously, the Alfven speed drops to zero at the null point, indicat- 
ing that the wave will never actually reach there, but the length scales (this can be 
thought of as the distance between the leading and trailing edges of the wave pulse) 
rapidly decrease, indicating that the current (and all other gradients) will increase. 
As an illustration, consider the wavefront as it propagates down the z— axis along 
X = y = 0. From Equation (22), the leading edge of the wave pulse is located at 
a position z = zoe"^'^^^^*, when the wave is initally at z = zq. If the trailing edge 
of the wave pulse leaves z = zq at t = ti then the location of the trailing edge of 
the wave pulse at a later time is Z2 = zoe"^*^^^^^*"*^'. Thus, the distance between 

the leading and trailing edges of the wave is Sz = ^oe"'-'^^^-'* ^e'*^^^-**^ — 1^ and this 

decreases with time, suggesting that all gradients will increase exponentially. 

We can also find analytical solutions for the velocity and polarisation of the fast 
wave. V • Bo = along x = y = implies Vz = 0, and hence using Equation (22) we 
obtain: 



(0,0, 



Using these forms in Equation (20) gives: 



J- = -(7fl)T(-"«'"-»)'* = -I7^'"""<--"-°"* ' 

where we have substituted for z from Equation (22). Thus, along the 2:— axis current 
builds up exponentially: [ji| ~ 2^^ ~ e'-"^''"^^*. Comparing to Equation (21) we see 
that this exponent is the same as that of our theoretical maximum current build up 
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Figure 7. (e = 1) Surfaces of constant at four values of t, showing the behaviour of the 
(initially planar) wavefront that starts at —0.25 < x < 0.25, —0.25 < y < 0, zg = 1 (blue) and 
—0.25 < a; < 0.25, < y < 0.25, = 1 (red). The (arbitrary) colouring has been added to 
aid the reader in tracking the wave behaviour. 

(under these initial conditions a = loq/ {e+ 1)). The coefficient is shglitly smaller 
than our theoretical maximum, but this is most likely because the limits we assumed 
for Equation (35) (see Appendix C) were not very strong. 

Secondly, for a planar fast wave starting at y = yo (Section 4.2), we can perform 
the same analysis along x = z = 0. Using the appropriate initial conditions (Equation 
19) and following the same analysis as above, we obtain: 

y = yoe-'\ v=iv,,0,v,)e''^° , k = (O, goe'*, O) 
=>Ji = e [vz,0,-vx)e^° . 

<^yo 

Hence, we have exponential current build up: [jij ~ y^^ ~ e'^*. Again, this exponen- 
tial build up is within our theoretical limits (Equation 21). 

5. Alfven Wave 



We now consider the second root to Equation (14) which corresponds to the Alfven 
wave. Hence, we assume uj'^ ^ |Bo| |k| and simplify Equation (14) to: 

J^{(t),x,y,z,p,q,r) = tj^ - (Bq • k)^ 

= - [xp + eyq - (e + 1) zrf = . (24) 
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Figure 8. (e = 1) Left: Ray paths for fluid elements that begin at points —0.25 < a;o < 0.25 
along j/o = 0, zq = 1 (indicated in black), and —0.25 < y < 0.25 along xq = 0, zq = 1 
(indicated in red) after a time t = -k 12. The ray path from xq = yo = is indicated in green 
and corresponds to the spine fieldline. Right: Projection of ray paths onto the xj/— plane (red 
indicates y > 0, blue y < and green y = 0. 



Charpit's Equations relevant to Equation (24) are: 



ds 







dt 
d^ 



dx 



, — = , — =eqE, 
ds ds ds 



dy 
ds 
dr 



-eye , ^ = {e+l)zi , 



(e + l)rC 



(25) 



where ^ = xp + eyq — {e + 1) zr. Thus, we can see that (j) = constant = and 
uj = constant = ujq. In addition, t = ujs, where we have set t = at s = 0. 



5.1. Planar Alfven Wave starting at = 1 



We now solve Equation (25) as before, subject to the initial conditions: 

(j)0 = , uo = 2tt , -l<a::o<l, -1<2/0<1, -20 = 1, 

po = , go = 0, ro=LJo/[{e + l)zo] , (26) 

where we have (arbitrarily) chosen ujq = 27r and (po = 0. This corresponds to a 
planar Alfven wave initially at z = zq. 

We can see the behaviour of the Alfven wavefront in Figure 7 (we have plotted 
surfaces of constant </> as in Section 4.1). We have also only plotted the wavefronts 
originating from —0.25 < xo,yo < 0.25 so as to better illustrate the wavefront 
evolution. We can clearly see that the initially planar wavefront expands (in the 
X?/— plane) as it approaches the null point, and keeps its original shape (i.e. planar 
and no rotation). The Alfven wave eventually accumulates along the fan plane, and 
never enters the z < domain. 

In Figure 8: Left, we can see the ray paths for fluid elements that begin at points 
-0.25 <xo< 0.25, yo = 0, zo = 1 and -0.25 < yo < 0.25, xq = Q, zq = 2, after a 
time t = tt/2. Here, we see that the fluid elements travel along and are confined to 
the fieldlines they start on, i.e. the Alfven wave spreads out following the fieldlines. 
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Figure 9. (e = 1) Surfaces of constant at four values of t, showing the behaviour of the 
(initially planar) wavefront that starts at —1 < a:o < 1, lyo = 1 and —0.5 < 20 < (blue) and 
— l<xo<l, 1/0 = 1 a-iid < 20 < 0.5 (red). The (arbitrary) colouring has been added to aid 
the reader in tracking the wave behaviour. The green dot indicates the position of the magnetic 
null point. We have imposed maximum and minimum values of unity in the 2— direction, purely 
for illustrative purposes. 

This explains the expansion of the wavefront seen in Figure 7. A similar effect was 
seen in the 2D case (McLaughlin and Hood, 2004). As noted for the wavefront, all 
the elements have travelled a different distance along their respective fieldlines but 
still form a planar wave. This is explained in section 6. 

5.2. Planar Alfven Wave starting at ?/o = 1 



We again solve Equation (24) but now subject to the initial conditions: 

(j)0 = , Uo = 2tt , -1 < XQ < 1 , yo = 1 , -1 < 20 < 1 , 

Po = , go = t^o/ (eyo) , ro = . (27) 

This corresponds to an Alfven wave being sent in from the side boundary (along 

y = yo)- 

We can see the behaviour of the Alfven wavefront in Figure 9 (surfaces of constant 
(/>). We have plotted the wavefronts starting at —1 < a;o < 1, yo = 1 and —0.5 < 
zo < 0.5 in order to more clearly show the Alfven wave propagation. We can see that 
the (initially rectangular) wavefront expands in the 2— direction but is also squeezed 
in the a;— direction as it approaches the null (i.e. as y decreases). We have imposed 
maximum and minimum values of unity in the z— direction, purely for illustrative 
purposes. The Alfven wave, and hence the wave energy, eventually accumulates along 
the spine. Again, the wave remains planar as it propagates. 

In Figure 10: Left, we can see the ray paths for fluid elements that begin at points 
— 1 < xo < 1, yo = 1 and at z = —0.25,0,0.25 after time t = 27r, where we have 
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■ ■ . ■■■■■ -1.0 -0.5 0.0 0.5 1 .0 

-Y 

Figure 10. (e = 1) Left: Ray paths for fluid elements that begin at points —1 < xq < 1, 
yo = 1 and zq = ±0.25 (indicated in red) and zq = indicated in black) after time t = 2it. We 
have imposed maximum and minimum values of unity in the 2— direction, purely for illustrative 
purposes. Right: Ray paths for fluid elements that begin at points —1 < zg < I in the yz— plane 
along X = (blue indicates starting points of —1 < zq < 1 in divisions of 0.1, red indicates 
20 = ±0.25, black indicates zq = 0). We have also plotted y — > —y to aid the comparison 
between the left and right figures. 



imposed maximum and minimum values of unity in the z— direction (again purely 
for illustrative purposes). The ray paths in the fan plane all focus towards the null, 
which is expected as they follow the fan-fieldlines. In contrast, the fluid elements on 
fieldlines above and below the fan plane propagate away from the null point, but are 
simply following their respectively fieldlines. This is also clearly seen in Figure 10: 
Right, which shows various ray paths in the yz— plane along a:: = 0. This behaviour 
explains the narrowing and stretching effect seen in Figure 9: the Alfven wave crosses 
the fan plane in this scenario and thus travels along the radially converging fan 
plane fieldlines. Meanwhile, the stretching effect comes from the diverging fieldlines 
the wave initially crosses. This work highlights the importance of understanding the 
magnetic topology of a system. 



6. Analytical Solution for Alfven wave 

We can also solve Charpit's Equations for the Alfven wave (25) analytically. Firstly, 
let us consider a planar wave starting at z = zq- Using the appropriate initial 
conditions (Equation 26), we find: 

df d , / 1 ^ X 

— = —{xp + eyq-{e + l)zr)=0 
as as 

^ = XQpo + eyoqo - (e + 1) zoro = -ujq (28) 

where ^ = xp + eyq — (e + 1) zr as before, and where the values of xq, po, yo, go, 
ZQ, ro and the sign of luq are taken from Equation (26). Thus, Equation (25) can be 
solved analytically: 

p = poe , q = QQe , r = roe' , 

X = xae , y = yoe , z = ZQe ' ' . (29) 
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where t = ujos. This sohition is vahd for all e. For e = 0, we recover the 2D solution 
of McLaughlin and Hood (2004). 

From these Equations we can see why an initially planar wave remains planar: if 
p and q are initially zero (po = 90 = 0) they remain zero for all time. In addition, 
z is independent of starting position xa and yo. Thus, after a given time, different 
elements have travelled different distances along their respective fieldlines, but all 
have the same z, i.e. all remain planar if originally planar. In addition, it can be 
shown that the volume occupied by the Alfven wave pulse is conserved (Appendix 
D). 

Consider a circular wavefront at z = zq, such that xq + i/q = r^, where r is some 
chosen radius. Let xi, yi, z\ represent the position of the wavefront after some time 
t. Thus, the change in length scales ((5x) can be represented as: 

(5a; = (xo - xi) e* , Sy = {yo - yi) e^* , 5z = {zq - zi) e^'-"'^^^* . 

Thus, the wave eventually accumulates along the fan plane, i.e. 5x oo, Sy — > oo, 
5z — > 0. Furthermore, the circular wavefront evolves as: 

2, 2 2 ^ f , f y V 2 

i.e. the circle becomes an ellipse (with semimajor-axis in the direction of x if < 
e<l,j/ife>l and remains circular for e = 1). Thus, the wave only accumulates 
over the whole fan plane for e = 1 and instead accumulates along a preferential axis 
for e ^ 1. 

Charpit's Equations (Equation 25) can also be solved using the initial condi- 
tions for a planar wave starting at y = yo, i-s. Equation (27). Following the same 
techniques above, we see that the length scales evolve as: 

5x = {xQ - xi) e~* , Sy = (yo - yi) e'"^ , 5z = {zq - zi) e^"'^^^^ . 

In this case, the wave eventually accumulates along the spine, i.e. 5x — > 0, (5y — > 0, 
5z — > oo, for all values of e. As before, an initially circular wavefront becomes 
elliptical for < e ^ 1, and evolves according to: 

2,2 2 ^ / a; \2 / 2 n2 2 



5(e + l)t 



6.1. Wavevector and Velocity 



From Section 3, we know that for the Alfven wave v ■ Bo = 0, v ■ k = and 
V • (Bo X k) 7^ 0. Consider a planar wave starting at 2 = zq, using Equation (29) we 
obtain: 

VxXoe +VyeyQe —Vz{e + l)zoe ^ ' = 

-t , -et , U+l)t n 

VxPoe + Vyq^e + VzVQe^ = 

Using the initial conditions from Equation (26), po = 90 = and so Vz = 0. Thus: 

k=(o,0,roe(^+i)*) , v = i;, (^-ge('-i)*,l,o)e'^° . (30) 



McLaughlin3DWKB.tex; 2/02/2008; 12:56; p. 18 



Application of WKB Theory to 3D MHD 



Thus, the angle between Vx and Vy changes with time. There is one special case: 
for e = 1, we have Vx/vy = —y/x = — tan^?. Recall in cylindrical coordinates 
Vx = Vr cos 6 — vg sin 9, Vy = Vr sin 6 + vg cos 9, and so we must have Vr = Q and 
vg 7^ 0. Hence, for e = 1 we have circular rotation of the fieldlines. 

Similarly, for a planar Alfven wave starting at y = yo (Equation 27) and using 
the same derivation as above, we obtain: 

k=(0,9oe^*,0) , v = ^;. (^(e + l)i^e(^+2)t^o^^j^i0o 
Finally, for a planar Alfven wave starting si x = xq, we obtain: 

k=(poe*,0,0) , v = .Jo,i + l^e(2^+i)*,l')e^*o . 

\ <^ yo J 

These velocity and polarisation solutions will be used in the next section. 



6.2. Current build up 



Recall from Section 2 that the perturbed electric current is given by ji = V x Bi. 
Now that we have an analytic solution for v we can solve Equation (8) for Bi . Hence, 
ji can be found: 

d 

— Bi = Vx(vxBo) ^ -tJoBi = k X (v X Bo) = (Bo • k)v - (k- v)Bo , 
^ ji = ik X Bi = -i (k X v) (Bq ■ k) /ujq , 

where we have made use of v • k = 0. Let us first consider a planar wave starting at 
z = ZQ. Using the forms of v and k from Equation (30) gives: 

ji = -i(k X ■v)i/uQ = - Vy [xQS^^ > ,yQe , e ^" , (31) 

(e + l)XQZQ \ / 

where ^ = Bo • k = — wo and luq = {e + 1) z^r^ from Equation (28). 

Thus, we have an exponential build up of jx and jy in our system. For e = 0, this 
reduces to jx ~ e* as found by McLaughlin and Hood (2004). We also see that the 
current build up is the fan-plane. 

Similarly, for a planar Alfven wave starting at y = yoi we obtain: 

j, = _Jfl^„ r,„e^*,0,-(e+l).oe^(^+i)*)e'*" , (32) 

where ^ = Bo ■ k = ojo and loq = eyoqo (from Equation 27). 
Finally, for a planar wave starting at a: = a;o, we obtain: 

Ji = 0, -eyoe ,(e + l)2oe^^ e ^" , (33) 

exoya \ I 

where ^ = Bo ■ k = wo and ix^o = a^opo- For e = 1, Equation (33) is the same as 
Equation (32) tmder the transformation (x,y) — > {—y,x). We can see that for both 
Equations (32) and (33), the current build up is predominately along the spine. 
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7. Conclusion 

We have demonstrated how the WKB approximation can be used to help solve the 
linearised MHD Equations. Using Charpit's Method and a Runge-Kutta numerical 
scheme, we have demonstrated this technique for a general 3D potential magnetic 
null point (parameter e). Under the assumptions of ideal and cold plasma, we have 
considered two types of wave propagation: fast magnetoacoustic and Alfvenic. 

For the fast magnetoacoustic wave, we find that the wave experiences a refraction 
effect towards the magnetic null point. The magnitude of the refraction is different 
for fluid elements approaching the null from various directions and is governed by the 
Alfven speed profile, v\ = x'^ + e^y^ + (e + 1)^ z'^ (in non-dimensionalised variables) 
and it is this different dependence on x, y and z that lead to different strength 
refraction effects. However, for all e the main result holds: the fast wave accumulates 
at the null point. 

In both Sections 4.1 and 4.2, the fast wave, and thus the wave energy, accumulates 
at the null point. The fast wave cannot cross the null because the Alfven speed there 
is zero. Thus, the length scales between the leading and trailing edges of wave pulses 
will decrease indicating that the current (and all other gradients) will increase. In 
Section 4.3, we calculated theoretical limits of the current build up and found that it 
was bounded by two exponentially growing functions. Moreover, it was shown that 
for a fast wave starting at z = za, |ji| r-j ~ e^*^^^-**, and for fast wave starting 
at y = j/o: ~ y~^ ~ e*^*. Hence, no matter how small the value of the resistivity 
is, if we include the dissipative term then eventually the 77V^Bi term in Equation 
(8) will become non-negligible and dissipation will become important. In addition, 
since ji grows exponentially in time, diffusion terms will become important in a time 
~log77; as found by Craig and Watson (1992) and Craig and McClymont (1993). This 
means that linear wave dissipation will be very efficient. Thus, we deduce that 3D 
null points will be the locations of wave energy deposition and preferential heating. 

We find that the Alfven wave propagates along the fieldlines, and that an Alfven 
wave fiuid element is confined to the fieldline it starts on. For the Alfven wave 
approaching the null point from above (planar wave starting at x = zq) the wave 
accumulates along the fan plane. For an Alfven wave approaching from the side 
(propagation initially perpendicular to the spine) the wave accumulates along the 
spine. This behaviour is in good agreement with the results of Pontin and Galsgaard 
(2007) and Pontin, Bhattacharjee, and Galsgaard (2007), but the method we present 
here clearly illustrate why this occurs, e.g. by following the ray paths in Section 5.2, 
it is clear why an Alfven wave generated crossing the fan plane must accumulate 
along the spine. 

Furthermore, we found an analytical solution for the Alfven wave. From this we 
were able to show that the Alfven wave rotates the fieldlines, the volume occupied by 
the wave pulse is conserved and that the associated currents build up exponentially. 
For a wave starting at z = zp, the currents build up along the fan plane, and jx 
and jy grow as e^'^'*'"'^-'* and e^*^*, respectively. Thus, resistive effects will eventually 
become non-negligible in a time ~ log 77 . For a wave starting at 2 = zq, the value 
of e determines where the preferential heating will occur: fan-plane (e = 1), along 
the X— axis (0 < e < 1) or along the y— axis (e > 1). In contrast, an Alfven wave 
starting at a:: = or ?; = yo will lead to preferential heating along the spine. 

All of the work described here highlights the importance of understanding the 
magnetic topology of a system, specifically the location of the spines and fans for 
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a 3D mill point. It is at these areas where preferential heating will occur, i.e. these 
areas are where the wave energy accumulates. In addition, it is of note that for both 
the fast and Alfven waves, current builds up exponentially and thus diffusion terms 
will become important in a time that depends on logr/. This is all in good agreement 
with the 2D work of McLaughlin & Hood (2004; 2005; 2006a). 

It is also useful to make an order of magnitude estimate for the quantities pre- 
sented here, in order to gain a better understanding of the physical conclusions. Let 
us consider oiu e = 1 system to have characteristic length L = 10 Mm, B = 10 G, 
po = 10~ kgm~ , fi = Anx 10~ Hm~ and 77 = 1 m s~ . This gives a characteristic 
speed of i; = 892km s~^, a characteristic time t = L/v = 11.2seconds, frequency 
LOO = 2n/t = 0.56 Hz, wavelength A = v/uq = 1.59 Mm and j = B / ij,L = 8x lO^^A. 
Thus, for the planar fast wave starting at zq = L and considering the behaviour along 
X = y = 0, we find that after a time t = I second, we have built up a current of 

0. 3mA (Equation 23). We can also estimate the time it takes for resistive effects to 
become important. We assume dB/dt « TyV^B uB f« r]B/{5z)'^, where 6z is the 
distance between the leading edge and trailing edge of our wave pulse, and we take 
5z = A, where the form of z is given by Equation (22). We find that resistive effects 
become non-negligible in a time t w log (ti;A^/ 77)/4 = 7 seconds. For comparison, 
after t = 7 seconds, our wave has built up a current of 0.87 mA and has travelled 
a distance of 7. 13 Mm. The Alfven wave is degenerate with the fast wave along the 
spine and so has the same estimates as above (under identical conditions). 

The 3D WKB technique described in this project can also be easily applied to 
other magnetic configurations, e.g. 3D dipole, and we hope that this paper has 
illustrated the potential of the technique. In addition, it is possible to extend the 
work by dropping the cold plasma assumption. This will lead to a third root of 
Equation (14) which will correspond to the behaviour of the slow magnetoacoustic 
wave. 

We conclude this paper with some caveats concerning the method presented here, 

1. e. if modellers wish to compare their work with a WKB approximation, it is es- 
sentialy to know the limitations of such a method. Firstly, in linear 3D MHD, we 
would expect a coupling between the fast and Alfven wave types due to the geometry. 
However, under the WKB approximation presented here, the wave sees the field as 
locally uniform and so there is no coupling between the wave types. To include the 
coupling, one needs to include the next terms in the approximation, i.e. the work 
presented here only deals with the first-order terms of the WKB approximation. 

Secondly, note that the work here is only strictly valid for high-frequency waves, 
since we took and hence a; = c^t to be a large parameter in the system. The 
extension to low frequency waves is considered in Weinberg (1962). 

Finally, the WKB approximation becomes degenerate at the points Vj^ = Cs, 
i.e. regions where the Alfven speed and sound speed are equal. Thus, the WKB 
method in the form presented here cannot be used to investigate mode conversion 
{e.g. see McLaughlin and Hood, 2006b) and, as mentioned above, the next terms 
in the approximation are needed. Alternatively, work is underway to overcome this 
degeneracy using the method developed by Cairns and Lashmore-Davies (1983) to 
match WKB solutions across the mode conversion layer (layer where Vj^ = Cs). The 
results of such work in ID can be found in McDougall and Hood (2007). 
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Appendix 

A. k parallel to Bq 

In this appendix, we address the scenario k = ABo in which the vectors of our three- 
dimensional coordinate system (Bo,k, Bq x k) are no longer linearly independent. 
To do this we consider the following Equation: 

^vi = c^V (V • vi) + {V X [V X (vi X Bo)]} X Bq , (34) 

which is derived in the same way as Equation (12) but without assuming a cold 
plasma. Under Cs = 0, Equation (34) reduces to Equation (12). Thus, assuming 
k = ABo and applying the WKB approximation (Equation 13) to Equation (34) 
gives: 

■p> 

w^v = (k • v) k + {k X [k X (v X Bo)]} x — 

MPo 

= (k • v) k + (k • Bo)' — - (k . Bo) (v • Bo) — 

MPO MPO 

- (k.Bo)(k.v)^ + (k-v)|Bot' — 
MPO MPo 

= c^A (Bo-v)Bo + A J '-Bo v-A ^ ^ (v ■ Bo) Bo 

MPO MPO 

_ A2mi£(Bo-v)Bo + A2(Bo-v)^Bo 

PPQ PPQ 

= c^a' (Bo ■ v) Bo + X^va |Bo|' v - X^va (v ■ Bo) Bo 

2 IB P 

where va = and we have explicitly included /i and po- Thus, for v parallel to 

Bo, i.e. V = aBo, we have: 

(jj'aBo = c^A^Q |Bo|' Bo + A^i;^ |Bo|' qBo — A^w^a |Bo|' Bq 

^ 2 2 I, |2 

=Y- LU = |k 

So the longitudinal oscillations (since v || Bq || k) propagate at the sound speed, 
i.e. this is the dispersion relation for slow waves. 

For V perpendicular to Bq (v • Bq = 0), i.e. transverse oscillations, we have: 

uj'^'v± = )?va |Bo|' Lo^ = Va |k|' 

This is the dispersion relation for a transverse and incompressional Alfven wave (i.e. 
k II Bo J- v). However, it is also the dispersion relation for the fast magnetoacoustic 
wave propagating in the direction of the magnetic field. Thus, we cannot distinguish 
between these two wave types in this specific scenario. 

It is also worth noting that even though the coordinate system we considered in 
Section 3 is not linearly independent when Bo || k, the result. Equation (14), still 
holds. Under the assumption k = ABo, Equation (14) simplifies to: 

^ {<l},x,y,z,p,q,r) = (^uj-vaI^^^ =0 
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Figure 11. (e = 1/2) Ray paths for fluid elements that begin at (xq, j/o, zq) = (A) (0, 1,1), {B) 
{— 1, 1, 0), (C) (—1, 1, 1) and (D) (—0.5, 1, 0.5). This figure is also available as an mpg animation 
showing many more starting points. Red indicates the 3D ray path and black indicates the xy, 
yz and xz projections of this ray path onto the respective planes. The blue dot indicates the 
position of the magnetic null point. 

So we have a double root and the solution is degenerate, i.e. it is impossible to 
distinguish the waves under these conditions (in agreement with the work above). 

B. (e = 1/2) Planar fast wave starting at yo = 1 and xq = 1 

The ray paths corresponding to a planar fast wave starting at yo = 1 and starting 
at xo = 1 in the e = 1/2 magnetic configuration can be found in Figures 11 and 
12. Recall that the fast wave cannot cross the null because the Alfven speed there is 
zero. 



C. Limits on fast wave current build up 

Define i?^ = x'^+y^+z'^ and assume < e < 1. RecaU |Bo|^ = a;^+e^2/^ + (e + 1)^ z^. 

2 ^ 2 ^ 2 

This can be bounded above by (e + 1) R , since 1 < (e + 1) and < (e + 1) . 
Hence, using a similar lower bound, we have: 

e^R^ < \Bq\^ < {e + lf R^ , (35) 
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Figure 12. (e = 1/2) Ray paths for fluid elements corresponding to a fast wave starting at 
xo = 1- We have made the transformation {x,y) {—y,x) in order to more easily compare 
with Figure 11. Thus, this figure shows ray paths that begin at {yo,xo,zo) = (A) (0,1,1), 
(B) (-1,1,0), (C) (-1,1,1) and (D) (-0.5,1,0.5). This flgure is also available as an mpg 
animation in the electronic edition of Solar Physics, showing many more starting points. Red 
indicates the 3D ray path and black indicates the xy, yz and xz projections of this ray path 
onto the respective planes. The blue dot indicates the position of the magnetic null point. 



These limits can also be understood physically: Recall that constant values of |Bo|^ 
defines an ellipsoid. Since we assume < e < 1, the largest distance from the centre 
to any edge of the ellipsoid is (e + 1) R, and the smallest distance is eR. Thus, 
physically we have encased our ellipsoid inside two spheres of radii eR and (e + 1) i?. 
From Equation (16) we have: 

^R' = 2x^+2y^ + 2z^ 
ds ds ds ds 

= -2 {xp + yq + zr) I'Bol'^ = -2 {xopo + yoqo + zoro) I'Bol'^ , 

where we have used the conserved quantity from Equation (17). Define a = xqpo + 
yoqo + zoro- Thus, from Equation (35) we have the inequality: 

-2a (e + 1)^ R^ < —R^ < -2ae^R^ , 
ds 

We can integrate and invert this inequality to obtain: 

Roe ^ ' < R < Roe =^ -^e - ^ - ^^^> 
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where Rq is a constant that depends upon starting position: Rq = xq + i/q + zq. 
Hence, inverting Equation (35) and combining it with Equation (36) gives: 

e < < 77 < -7[ 7T < — :: ^e ^ ^ ' 



Finally, we recall t = ujqs and thus 



1 Jat^t/iJo ^ 1 ^ _J_ 2a{e+iyH/uJo 



{e + lfRl -\Bof-e^Rl 



D. Volume 

Assume we generate an initially rectangular wave pulse of volume Vb = {xi — X2) x 
{yi ~ 2/2) X {zi — Z2), where xi, X2, yi, y2, zi and Z2 define the starting points at 
the edges of our domain. The wave will evolve according to Equation (29) and thus, 
after travelling distance s along the characteristic curve, will occupy a volume: 

end = yxie - X2e ) X (yie - y2e ) x (216 ^ ' - Z2e ^ ' j 
= {xi - X2) X {yi - j/2) X {zi - 22) = Vb 

Thus, volume is conserved for an Alfven wave in this system. 
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